Interdisciplinary 

Mathematics 

Institute 

2014:08 

A  3D  Hydrodynamic  Model  for 
Cytokinesis  of  Eukaryotic  Cells 


Jia  Zhao  and  Qi  Wang 


College  of  Arts  and  Sciences 
University  of  South  Carolina 


Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
VA  22202-4302  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  a  penalty  for  failing  to  comply  with  a  collection  of  information  if  it 
does  not  display  a  currently  valid  OMB  control  number 

1.  REPORT  DATE 

2014  2.  REPORT  TYPE 

4.  TITLE  AND  SUBTITLE 

A  3D  Hydrodynamic  Model  for  Cytokinesis  of  Eukaryotic 


6.  AUTHOR(S) 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

University  of  South  Carolina, Interdisciplinary  Mathematics 
Institute, Columbia, SC, 29208 

9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 


12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  unlimited 

13.  SUPPLEMENTARY  NOTES 

14.  ABSTRACT 

In  the  late  stage  of  the  mitotic  cycle  of  eukaryotic  cells  ,  cytokinesis  ensues  where  a  parent  cell  replicates  its 
nucleus  with  the  necessary  genetical  substances  (i.e.,  DNAs  and  chromosomes)  and  splits  into  two  similar 
offspring  cells.  This  mitotic  process  involves  complex  chemical,  biophysical  and  mechanical  processes 
whose  details  are  just  begin  to  be  unfolded  experimentally.  In  this  paper,  we  propose  a  full  3-D 
hydrodynamical  model  using  a  phase  eld  approach  to  study  the  cellular  morphological  change  during 
cytokinesis.  In  this  model,  the  force  along  the  contracting  ring  or  cytokinetic  ring  induced  by  remodeling  of 
actin-myosin  lament  on  cell  cortex  layer  at  a  division  plane  of  the  parent  cell  during  cytokinesis,  is 
approximated  using  a  proxy  force  anchored  on  the  newly  formed  nuclei.  The  symmetric  or  asymmetric  cell 
division,  i.e.  a  parent  cell  dividing  its  cytoplasm  to  produce  two  equal  sized  or  unequal  sized  offspring  cells, 
is  simulated  numerically  with  the  model.  Our  numerical  results  show  that  the  location  of  the  division  plane 
and  the  contracting  force  along  the  cytokinetic  ring  on  the  division  plane  is  essential  for  the  cell  division.  In 
addition,  our  numerical  study  also  shows  that,  during  cytokinesis,  surface  tension  of  the  cell  membrane 
also  contributes  to  this  process  by  retaining  the  morphological  integrity  of  the  offspring  cells.  This  model 
and  the  accompanying  numerical  simulation  tool  provides  a  solid  framework  to  build  upon  with  more 
sophisticated  whole  cell  models  to  probe  the  cell  mitotic  process. 

15.  SUBJECT  TERMS 


16.  SECURITY  CLASSIFICATION  OF: 

17.  LIMITATION  OF 

18.  NUMBER 

19a.  NAME  OF 

ABSTRACT 

OF  PAGES 

RESPONSIBLE  PERSON 

a  REPORT 

unclassified 

b  ABSTRACT 

unclassified 

c  THIS  PAGE 

unclassified 

Same  as 
Report  (SAR) 

21 

3.  DATES  COVERED 

00-00-2014  to  00-00-2014 

5a.  CONTRACT  NUMBER 

5b.  GRANT  NUMBER 

5c.  PROGRAM  ELEMENT  NUMBER 

5d.  PROJECT  NUMBER 

5e.  TASK  NUMBER 

5f.  WORK  UNIT  NUMBER 

8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 

10.  SPONSOR/MONITOR'S  ACRONYM(S) 

11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 


Standard  Form  298  (Rev.  8-98) 

Prescribed  by  ANSI  Std  Z39-18 


A  3D  Hydrodynamic  Model  for  Cytokinesis  of 

Eukaryotic  Cells 

.Jia  Zhao  *and  Qi  Wang  t 

Abstract 

In  the  late  stage  of  the  mitotic  cycle  of  eukaryotic  cells  ,  cytokinesis  ensues  where  a 
parent  cell  replicates  its  nucleus  with  the  necessary  genetical  substances  (i.e.,  DNAs  and 
chromosomes)  and  splits  into  two  similar  offspring  cells.  This  mitotic  process  involves 
complex  chemical,  biophysical  and  mechanical  processes  whose  details  are  just  begin 
to  be  unfolded  experimentally.  In  this  paper,  we  propose  a  full  3-D  hydrodynamical 
model  using  a  phase  held  approach  to  study  the  cellular  morphological  change  during 
cytokinesis.  In  this  model,  the  force  along  the  contracting  ring  or  cytokinetic  ring, 
induced  by  remodeling  of  actin-myosin  filament  on  cell  cortex  layer  at  a  division  plane 
of  the  parent  cell  during  cytokinesis,  is  approximated  using  a  proxy  force  anchored  on 
the  newly  formed  nuclei.  The  symmetric  or  asymmetric  cell  division,  i.e.  a  parent  cell 
dividing  its  cytoplasm  to  produce  two  equal  sized  or  unequal  sized  offspring  cells,  is 
simulated  numerically  with  the  model.  Our  numerical  results  show  that  the  location  of 
the  division  plane  and  the  contracting  force  along  the  cytokinetic  ring  on  the  division 
plane  is  essential  for  the  cell  division.  In  addition,  our  numerical  study  also  shows 
that,  during  cytokinesis,  surface  tension  of  the  cell  membrane  also  contributes  to  this 
process  by  retaining  the  morphological  integrity  of  the  offspring  cells.  This  model  and 
the  accompanying  numerical  simulation  tool  provides  a  solid  framework  to  build  upon 
with  more  sophisticated  whole  cell  models  to  probe  the  cell  mitotic  process. 


1  Introduction 

Cell  is  the  fundamental  unit  in  all  living  organisms  since  animals  and  plants  are  all  made 
up  of  cells  of  different  varieties.  The  study  of  cells  is  therefore  an  essential  part  of  research 
in  biological  science  and  medicine.  Among  many  functions  of  a  cell,  one  important  function 
is  the  cell’s  ability  to  reproduce,  which  is  also  known  as  cell  proliferation.  Because  of  cell 
proliferation,  all  living  organisms  can  grow  and  acquire  their  life.  Given  its  unique  role  played 
in  biological  systems,  cell  study  has  been  the  focal  point  of  biology  for  centuries.  With  the 
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advancement  of  experimental  technologies  today,  more  cell  functions  and  micro-structures 
that  endow  the  functions  have  been  uncovered,  which  reveal  an  amazingly  complex  cosmo 
on  its  own. 

One  fundamental  part  in  cell  study  is  the  cell  reproductive  cycle,  where  a  parent  cell 
undergoes  a  sequence  of  intercellular  transformations  and  eventually  divides  into  two  or 
more  offspring  cells.  For  prokaryotic  cells,  the  cell  proliferation  process  is  called  binary 
division  or  binary  fission,  which  is  the  primary  method  for  reproduction.  For  eukaryotic 
cells,  it’s  called  cell  mitotic  process  or  mitosis.  At  the  late  stage  of  the  cell  mitotic  process 
for  eukaryotic  cells,  after  the  nucleus  has  been  divided  and  chromosomes  separated,  the  cell 
division  process  is  also  called  cytokinesis.  For  eukaryotic  cells,  cell  division  is  a  much  more 
complicated  process  than  the  division  of  prokaryotic  cells.  Despite  of  extensive  studies  on 
cells  over  the  years,  lots  of  important  mechanisms  in  living  cells  still  remain  mysteries.  A 
lot  more  need  to  be  discovered  experimentally  as  well  as  theoretically  in  order  for  one  to 
fully  understand  how  a  tiny  cell  really  works,  ft’s  commonly  agreed  upon  that  cell  division 
is  not  simply  driven  by  a  single  mechanism,  but  rather  many  factors  which  are  intertwined. 
As  one  of  the  most  spectacular  part  of  the  cell  cycle,  any  single  step  goes  wrong  may  lead 
to  a  catastrophe  or  failure,  which  may  lead  to  an  unwelcome  outcome  for  instance  cancer. 
Thus,  a  detailed  understanding  on  cytokinesis  can  be  of  great  benefits  for  understanding 
many  diseases.  For  more  details,  readers  are  referred  to  the  insightful  review  article  on  this 
topic  in  [E3].  Readers  can  also  find  comprehensive  review  materials  for  animal  cytokinesis 
in  [□,  03]  and  for  cytokinesis  of  bacteria  in  [El].  Besides,  some  works  related  to  mechanical 
properties  of  cells  during  cytokinesis  such  as  material  properties  of  cells  and  sources  of  stresses 
can  be  found  in  [EH,  E0].  The  study  published  in  [0]  discusses  the  molecular  requirements 
for  cytokinesis  and  the  work  in  [El]  addresses  some  recent  advances  in  the  mechanism  of 
cytokinesis  in  animal,  yeast  and  plant  cells. 

Experimental  observations  have  provide  us  with  a  basic  picture  of  cell  mitosis.  For  eu¬ 
karyotic  cells,  at  the  beginning  of  the  cell  mitotic  process,  the  parent  cell  first  duplicates  its 
genetic  substances  and  then  forms  a  mitotic  spindle  consisting  of  microtubules  [EH].  Through 
a  cascade  of  signaling  processes  [E3],  the  actin  and  myosin  molecules  would  undergo  a  self- 
assembly  process  to  remodel  the  cell  cortical  layer,  a  layer  rich  in  actin-hlaments  located 
immediately  adjacent  to  the  cell  membrane  [ESI].  In  sync  with  the  elongation  of  the  mitotic 
spindle,  more  actin  and  myosin  molecules  ascend  to  a  ring  like  region  in  a  plane  roughly 
orthogonal  to  the  axis  of  the  mitotic  spindle  to  form  the  cytokinetic  ring  or  contractile  ring. 
The  plane  is  called  the  cleavage  plane  or  division  plane  [0].  As  more  actomyosin  molecules 
are  accumulated  along  the  cytokinetic  ring,  a  contracting  force  is  generated  directed  inward 


toward  the  axis  of  the  spindle  [EH]].  The  contracting  force  pushes  the  membrane  inward 


to  create  what  is  known  as  the  cleavage  furrow  [ESI] .  The  localized  activation  of  the  small 
GTPase  Rho  family  of  proteins  at  the  cell  division  plane  controls  the  position  of  the  contrac¬ 
tile  ring  [ESI].  When  Rho  is  specifically  activated  at  the  division  plane  within  the  cortex,  it 
promotes  actin  polymerization  and  myosin-2  activation  via  Rho  effector  proteins.  Rho-GTP 
promotes  actin  filament  assembly  and  myosin- 2  assembly  [0],  The  contractile  ring  is  a  dy¬ 
namic  structure,  in  which  F-actin  and  myosin-2  are  continuously  assembled  and  disassembled 
to  maintain  a  contracting  force  to  squeeze  the  cell  along  the  contractible  ring.  For  eukaryotic 
cells,  the  positioning  of  contractile  ring  and  cleavage  furrow  is  regulated  by  distribution  of 
long  astral  microtubles  (MT)  [El].  Besides,  it  has  been  shown  that  the  cell  geometry  can 
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influence  the  position  of  the  cell  division  plane  [tE3]  and  ultimately  affect  the  morphology  of 
the  offspring  cells. 

In  addition  to  the  vast  amount  of  data  from  experimental  observations,  theorists  have 
been  trying  to  come  up  with  testable  hypothesis  and  models  to  decipher  underlying  mech¬ 
anisms  that  control  the  cell  division  process.  In  [□],  the  author  proposed  a  mathematical 
model  for  cell  cleavage  for  the  sea  urchin  by  considering  chemotactic  motion  of  the  centro- 
somes.  In  [E3],  the  author  developed  a  mechanical  cell  division  model  by  a  level-set  approach, 
highlighting  the  furrow  thinning  trajectories.  In  [1213],  the  author  formulated  an  immersed 
boundary  approach  for  modeling  cell  division  focusing  primarily  on  the  mechanical  aspect  of 
cell  division,  which  has  been  extended  by  [E3]  to  study  a  single  axisymmetric  cell  growth  and 
division.  The  correlation  between  cell  shape  elongation  and  the  orientation  of  the  division 
axis  has  been  studied  in  [ESI].  Besides,  a  model  based  on  polymer  dynamics  to  study  spindle 
dynamics  during  cell  division  has  been  developed  in  [BO]. 

Given  the  theoretical  attempts  made  so  far,  there  has  been  very  little  theoretical  work 
done  to  resolve  cytokinesis  in  full  3D  using  hydrodynamic  models.  Hence,  in  this  paper, 
we  set  out  to  develop  a  3D  hydrodynamic  model  for  cell  division  systematically  by  a  phase 
field  approach.  Phase  field  methods  have  gained  tremendous  attention  in  last  few  decades, 
in  particular,  for  the  study  of  interfacial  problems.  The  advantage  of  phase  field  methods, 
compared  to  the  traditional  sharp  interface  method,  is  that  the  interface  is  naturally  incor¬ 
porated  in  the  model  formation  such  that  it  does  not  need  to  be  tracked  separately.  This 
reduces  the  computational  cost  tremendously.  In  the  phase  field  formation,  there  is  a  param¬ 
eter  controling  the  thickness  of  the  interface.  The  phase  field  model  converges  to  its  sharp 
interface  limit  [0]  as  the  thickness  of  the  interface  approaches  zero.  For  cells,  the  coarse¬ 
grained  membrane  and  other  interfaces  within  a  cell  do  not  need  to  be  very  sharp,  which 
makes  the  phase  field  formulation  a  good  choice  for  describing  the  multiphasic  structure  in 
a  cell. 

Our  objective  in  modeling  cell  division  is  to  by  develop  a  mathematical  model  systemat¬ 
ically  based  on  a  phase  field  approach,  where  the  cytoplasm,  cell  nucleus,  and  extra-cellular 
matrix  are  treated  as  three  different  multiphasic  components  of  a  complex  fluid  mixture. 
Given  the  complexity  of  a  cell,  a  viable  model  must  be  built  step  by  step  by  incorporating 
more  complex  mechanochemical  details  incrementally.  The  first  important  step  is  therefore 
to  build  a  framework  to  simulate  the  morphological  change  of  the  cell  membrane  during 
cell  division  by  highlighting  the  most  dominant  mechanism  during  the  process.  The  most 
dominating  factor  that  one  has  identified  is  the  existence  of  the  contractile  ring  within  the 
cortex,  which  is  primarily  consisted  of  F-actins  and  myosins.  Modeling  the  detail  of  the 
cortex  remains  our  goal  for  this  project.  However,  the  effective  force  that  squeezes  the  cell 
during  cytokinesis  is  the  contractile  force  exerted  by  the  contractile  ring.  To  highlight  this 
mechanism  in  a  simplistic  way,  we  employ  a  proxy  force  in  lieu  of  the  contractile  force  gen¬ 
erated  by  the  actomysin  dynamics  in  the  cortex  in  the  model.  We  remain  committed  to 
develop  the  more  detailed  model  for  the  cortical  layer  and  the  contractile  ring  within  it,  but 
defer  it  to  a  future  paper.  The  governing  system  of  equations  in  the  model  consists  of  the 
transport  equation  for  each  phase  variable  in  the  form  of  the  volume  fraction  for  that  phase, 
continuity  equation  and  the  momentum  balance  equation.  We  discretize  the  governing  par¬ 
tial  differential  equations  and  solve  it  using  a  finite  difference  method  on  a  GPU  in  full  3D 
space  and  time.  Numerical  simulations,  which  will  be  discussed  below,  show  several  inspiring 
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morphologies  during  cell  division  that  match  well  with  the  experimental  observations. 

The  rest  of  the  paper  is  organized  into  three  sections.  In  section  two,  we  present  the 
multi-phase  hydrodynamic  model  for  the  whole  cell  in  a  simplified  picture.  It  is  followed 
by  a  brief  discuss  on  its  numerical  discretization  and  implementation  issues.  Finally,  we 
present  three  numerical  results  and  discuss  their  comparison  with  available  experimental 
observations. 


2  Mathematical  Formulations  of  MutliPhase- field  Mod¬ 
els 

In  this  section,  we  propose  a  mathematical  model  for  the  cell  division  process  given 
in  a  phase  field  formulation.  In  this  model,  we  treat  the  cell  and  the  surrounding  liquid 
environment  as  a  fluid  mixture  in  the  form  of  a  single  fluid  with  multiple  components.  The 
volume  fractions  of  cytoplasm,  nucleus  and  the  buffer  fluid  outside  the  cell  are  denoted  by 
0i,02  and  03,  respectively.  Here,  we  assume  the  buffer  is  a  viscous  fluid,  so  is  the  nucleus. 
The  cytoplasm  can  be  treated  as  a  viscoelastic  fluid.  In  this  paper,  we  also  treat  it  as  a  viscous 
fluid  within  the  time  scale  of  interest  for  simplicity.  rjt ,  i  —  1,2,3  are  used  to  represent  the 
respective  viscosities  for  the  three  phases,  respectively.  The  cell  membrane  (together  with 
the  cortical  layer)  is  the  level  set  defined  by  {0i  =  4  =  03}  and  the  membrane  of  the  nucleus 
is  the  level  set  defined  by  {02  =  \  =  0i}-  Notice  that  this  is  a  globally  multi-phase  while 
locally  binary  system  since  there  is  no  contact  between  nucleus  and  buffer  at  any  time.  The 
volume-average  velocity  and  density  for  this  fluid  mixture  is  defined  as 

3  3 

v  =  ^Vi’  P  =  ^2  l1) 

i= 1  i=  1 

where  v,,  /y  is  the  effective  velocity  and  density  for  component  i,  i  —  1,  2,  3.  For  incompress¬ 
ible  materials,  we  enforce 

01  +  02  +  03  =  1-  (2) 

2.1  Thermodynamic  free  energy 

We  denote  the  domain  in  which  the  cell  resides  together  with  the  buffer  fluid  as  fh  The 
free  energy  of  this  mixture  system  is  proposed  as  follows, 

F=  f  fd x,  (3) 

Jn 

where  /  is  the  free  energy  density  function.  There  are  different  choices  for  the  free  energy 
density  function  of  the  three  phase  fluid.  Here  we  adopt  a  simple  one: 

/  =  i(7A||V01||2  +  7B||V02||2  +  7c||V03||2 

V  \  (4) 
+7l0{0|  +  720203  +  73010?  +  7l230?0|0lJ  • 
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An  alternative  is  given  by 


/  =  H  7a||  V0;l||2  +  7^11  V02||2  +  7c||  V03|p 


+7l0l(1  -  0l)2  +  7201(1  -  02)2  +  7301(1  -  03)2  +  7123010202 


(5) 


3  h 


where  7^,  7 7c  govern  the  strength  of  the  conformational  entropy  between  different  compo¬ 
nents,  71,72,73  control  the  strength  of  the  bulk/mixing  energy  for  each  pair  of  components, 
and  7123  is  a  Lagrangian  multiplier  to  penalize  the  coexistence  of  the  distinct  phases.  Our 
numerical  studies  show  that  these  two  choices  of  free  energy  density  functions  yield  qualita¬ 
tively  the  same  results.  Thus,  in  our  study  presented  in  this  paper,  we  choose  the  second, 
i.e.  equation  (0).  We  remark  that  more  features  can  be  added  to  the  model  by  augmenting 
the  corresponding  free  energy. 

2.2  Transport  equations  for  biomass 

Given  the  specific  form  of  the  free  energy  density  (0),  we  assume  that  each  component 
in  the  fluid  mixture  is  convected  by  the  volume-averaged  velocity  as  well  as  transported  via 
the  osmotic  pressure.  Then,  the  transport  equation  for  phase  variables  are  given  as  follows 


dt<h  +  V  ■  (v^,)  =  V  ■  (  aij Vj—  )  +  g,,  i  =  1, 2, 3, 


(6) 


where  (cqj)  is  the  motility  matrix  and  Qi  are  the  reactive  terms,  respectively. 

Note  that  the  volume  fractions  add  up  to  1:  Yli=i  (pi  =  1  for  incompressible  mixtures. 
This  along  with  the  Onsager  reciprocal  principle  implies  [EES,  Id] 


^  ^  O/J  0,  Oiij  Oiji  • 

3= 1 

The  off-diagonal  mobility  coefficients  can  be  obtained  from  the  diagonal  coefficients 


(7) 


«12  —  ~  (®33  —  Ol\  1  —  0122)1  «13  —  ~{0l22  ~  «11  —  073),  OI23  ~  ^(“i1  _  tt22  _  a33),  (8) 


where  we  assume  the  diagonal  motility  parameters  have  the  following  form 

Ola  -^i0i  ( 1  0i), 


(9) 


with  A i  representing  the  strength  of  motility  parameters  for  the  corresponding  component  i, 
i— 1,2,3. 

By  incompressibility  condition,  it  follows  that 


1N  =  0- 


(10) 


i—  1 


At  the  onset  of  cell  division,  the  parent  cell  normally  increases  in  cytoplasmic  and  organelle 
volume  (the  G1  phase)  as  well  as  increases  in  genetic  materials  (the  G2  phase)  right  before 
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the  replication  during  the  S  phase.  This  process  can  be  modeled  by  proposing  the  reactive 
kinetics  for  the  time  rate  of  change  in  the  volume  fractions: 


g  1  =  Ci-f/l(£)0i03  -  C2-f/20)0102,  (H) 

$2  =  02^2(^)0102,  (12) 

03  =  -Cihfi(f)0i03.  (13) 


Here  Hi(t),i  =  1,2,  are  Heviside  functions  defined  by 


Hi(t) 


1,  t  t  1,  TJ  /,\  f  1,  t  ^2, 

0,  t  >  ti,  21  J  “  \  0,  t  >  t2, 


(14) 


where  H  is  the  critical  checkpoint,  at  which  the  parent  cell  has  just  duplicated  its  volume, 
f2  is  the  critical  checkpoint,  at  which  the  nucleus  has  doubled  its  volume.  Thus,  before 
fi  the  parent  cell  keeps  reproducing  its  cytoplasm  to  expand  volume.  After  it  doubles  its 
volume,  reproduction  of  cytoplasm  ceases  and  cell  proliferation  begins.  The  growth  kinetics 
is  assumed  to  be  originated  from  the  interface  between  the  buffer  and  cytoplasm  as  well  as 
the  nucleus  and  the  cytoplasm. 

We  remark  that  this  multi-phase  held  model  is  consistent,  in  term  of  reducing  to  the 
phase  held  model  with  less  physical  phases.  First  of  all,  when  the  component  i  is  not  present 
in  the  mixture  at  the  initial  time,  the  component  will  be  absent  from  the  system  during  the 
time  evolution  of  the  system,  i.e.  given  i  e  N, 


000)  =  0  — >  00i)  =  0,  Vf  >  0. 


(15) 


Secondly,  when  there  are  only  n  <  N  phases  in  the  model,  the  N-phase  model  naturally 
reduces  to  the  n-phase  model. 

Although  our  formulation  in  equation  (0)  is  consistent,  in  real  numerical  simulations,  we 
don’t  need  to  calculate  every  single  phase  due  to  the  incompressibility  constraint.  In  our 
later  discussion,  we  only  keep  track  of  the  transport  of  02  and  03,  since  0 1  can  be  obtained 
through  the  incompressibility  condition  (□). 


2.3  Continuity  and  momentum  equation 

In  order  to  close  this  system,  we  need  to  supplement  the  system  with  the  continuity 
equation  and  momentum  equation  for  the  fluid  mixture.  By  assuming  the  average  velocity 
in  the  fluid  mixture  solenoidal,  we  have 

p(0v  +  v  ■  Vv)  =  — Vp  +  V  •  r  +  Fe,  (16) 

Vv  =  0,  (17) 

where  p  is  the  hydrostatic  pressure  and  r  is  the  viscoelastic  stress  tensor  and  Fe  is  the  elastic 
body  force  yielding  the  surface  tension  of  the  cell  membrane  and  the  contractile  force  due  to 
actin-myosin  filaments.  The  viscoelastic  stress  tensor  is  proposed  as  follows, 

3 

r  =  2y>Di  +  r„  (18) 

i=l 
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where  the  first  term  sums  up  the  viscous  stress  from  each  component,  and  the  second  term 
re  is  the  elastic  stress  for  the  cytoplasm.  A  constitutive  equation  could  be  proposed  for  re 
relating  to  the  specific  structure  of  the  cortex  layer.  In  this  paper,  however,  we  set  re  =  0 
for  the  time  scale  of  our  interest  and  simplicity. 

For  the  elastic  body  force  or  the  interfacial  force  Fe,  we  adopt  the  surface  force  yielded 
by  the  variation  of  the  free  energy.  In  addition,  we  propose  a  proxy  force  mimicing  the 
cytokinetic  ring,  a  contractible  ring  responsible  for  the  cytokinesis  of  the  animal  cell,  which 
is  similar  to  the  one  used  in  [13 


F 


^  TT~ ^ +  ^i/dl0lV0i  +  Kfo ,/'<i202 V02 , 

j= 1  0CPj 


(19) 


where  the  first  term  is  the  interfacial  force,  due  to  the  material  change  at  the  interface 
beween  each  component  and  the  second  term  is  the  proxy  force  mimicing  the  normal  force 
generated  by  the  cytokinetic  ring,  consisting  of  actomyoisn  networks  immediate  within  the 
cell  membrane.  We  note,  this  proxy  force  provides  the  contractible  force  necessary  for 
dividing  the  cell.  As  we  alluded  to  in  the  introduction,  we  will  leave  the  model  of  tracking 
the  spatial  distribution  of  actin-hlament  networks  regulated  by  myosins  in  the  cytokinetic 
ring  for  a  future  work.  Specifically,  we  use  k  to  denote  the  mean  curvature, 

V0  ,  s 

K  =  -V.„,  n  =  -]wf.  (20) 

In  phase  field  formulation,  it  is  approximated  by 

K=j2Wv20-?0(l-0)  (1-20)1,  (21) 


where  £  is  the  thickness  of  the  interface,  and  fdi ,  i=l,2  denote  the  strength  of  the  cytokinetic 
ring  force,  which  by  following  the  idea  in  [03]  are  proposed  as  follows 


fdi  —  Hdi((p2 ,  t)- 


Idi 


X  -  cd  -  X  -  c2  +  £d 


fd2  —  Hd 2(03;  t)-. 


Id2 


X  -  Ci  -  X  -  C2  +  £d 


(22) 


where 


Hdl{4>2,  t) 


0, 

1, 


otherwise, 

02  =  0  >  t3. 


Hd2(4>3,t) 


0,  otherwise,  ,  s 

1,  03  =  0  &  t2  <  t  <  t3. 


are  heaviside  functions  to  restrict  the  force  on  cell  membrane  and  the  nucleus  membrane  at 
different  stages  of  the  cell  division  process,  respectively,  and  ci  and  C2  are  the  mass  centers 
of  the  separated  nuclei,  respectively,  i.e. , 

_  fQl  02(x)xdx  _  fna  02(x)xdx 
L,  02(x)dx  ’  Q  /n2  02(x)dx 

fix  and  fl2  are  the  domains  occupied  by  the  two  nuclei,  7 di,i  =  1,2  are  two  parameters 
characterizing  the  strength  of  fdi,  respectively,  and  ed  is  a  small  number  employed  here  to 


7 


avoid  the  singularity  of  fcjt .  Here,  t2  is  the  checkpoint  for  the  nucleus  to  begin  separating 
after  the  size  of  the  cell  has  been  doubled  and  £3  is  the  critical  checkpoint  when  the  nucleus 
has  been  divided  into  two  separate  nuclei. 

Notice  that  this  cytokinetic  ring  force  counters  the  force  clue  to  the  one  of  surface  tension 
along  the  cytokinetic  ring.  We  did  not  incorporate  this  cytokinetic  ring  force  into  the  ther¬ 
modynamic  free  energy  since  this  is  a  proxy  for  the  active  force  generated  by  actin-myosin 
filament  on  the  cytokinetic  ring  through  the  released  of  the  hydrolyzed  ATP.  In  the  current 
context,  it  is  not  a  potential  force. 


2.4  Dimensionless  governing  equations  in  the  three  phase  model 

We  denote  the  reference  time  scale  as  to,  reference  length  scale  as  h  and  reference  mass 
density  as  po-  Then,  the  variables  and  parameters  are  nondimensionalizsed  as  follows: 


fa  =  =  A,B,C,fs  =  =  1,2,3,123,  ^  =  ^,*  =  1,2,3 


ph4 


7 skTtl 


(24) 


p0h2 


riit0 


The  governing  equations  in  dimensionless  form  are  summarized  below. 


f  p(dtv  +  V  ■  Vv)  = 


Vv  = 


v  ■  (Ell  ffc&(Vv  +  VvT))  -  Vp  +  Fe, 

0, 


l  dt(j)i  +  V  ■  (v<j>i)  =  V  •  (E'Li  oiij  V f£))  +gh  i  =  1, 2, 3, 


where 


9i  =  -  c2H2(t)<j)1<f)2, 


(25) 


g2  —  c2H2(t)(f)i(f)2, 
g3  =  C\H\  (t)0i03, 


(26) 


f.  =  Eh.fvpj  +  E?.i^./«*VA. 

This  partial  differential  equation  system  will  be  discretized  using  a  finite  difference  method. 


3  Numerical  methods 

For  these  coupled  PDEs,  we  use  the  extrapolation  strategies  to  decouple  them.  That 
is,  for  each  time  step,  we  solve  the  momentum  equation  firstly  using  the  extrapolated  data 
for  the  biomass  and  functional  components.  Then  with  the  updated  velocity,  we  solve  the 
transport  equations  for  the  volume  fractions  and  the  functional  components  one  by  one. 
We  note  that,  in  the  following  context,  any  variable  with  an  overline  (•)"+1  represents  a 
second-order  extrapolation  from  n  and  n  —  1  steps. 


For  the  momentum  equation,  we  use  the  modified  Gauge-Uzawa  method  [010]  to  calculate 
the  average  velocity.  Recall  that  the  momentum  equation  is  given  by 


<9v 

p(—  +  v  ■  Vv)  =  V  ■  r  -  Vp  +  Fe. 


(27) 


By  adding  a  second  order  term  —  -U—  V2v  on  both  sides,  where  Rea  is  the  averaged  Reynolds 


number,  we  rewrite  the  momentum  equation  into 


,<9v  1 

p(w  +  v.vv)-- 


-V2v  =  —  Vp  +  Fe 


Re 


-V2v. 


(28) 


The  Gauge-Uzawa  method  is  given  by  three  steps  listed  below. 
1.  Prediction: 


P 


26t. 


-]+P" 


Vvn+1  +  |(V  •  (pn 


1))v' 


+  77-Vsn  +  Vpn  -  ^V2un+1  =  Rn+1  -  U— V2vn+1_e, 

ries  ^ea  “ea 

U”+1  •  Yi\y-Q  !  =  0, 


(29) 


2.  Projection: 


o"+! 


V++1)  =  V  •  u 


n+ 1 


0pn+1 I 

an  E=°>-U 


o, 


(30) 


3.  Correction: 


vn+1  = 


sn+1  = 


Un+1  + 


n+l 


sn  —  V  ■  u 


n+l 


^.n+1  „n  3-i/>n+1  ,  1  „n+l 

P  P  2 Si  Ren  ’ 


(31) 


where 


R"+1  = 

Z=1 


V  p"+1  V^”+1+5]  VVV1.  (32) 


n+l  —n+l 


z=l 


Z— 1 


and  e  =  0.05.  The  averaged  Reynolds  number  i?ea  is  computed  by 


min( 


1  1 

D  ’  77>  ’ 

itei  ite2 


(33) 


Where  i?,ei ,  i?e2  and  i?e3  are  the  Reynolds  numbers  for  cytoplasm,  nucleus  and  ECM  respec¬ 
tively.  Here  s°  =  0  and  v1,  s1,  <+  c 1  are  computed  by  a  first  order  scheme. 

For  the  multiphase  field  model,  the  coupled  Cahn-Hilliard  equations  could  be  solved 
simultaneously,  using  preconditioned  bicgstab  strategy.  However,  we  need  to  extrapolation 
the  nonlinear  terms,  which  results  in  a  coupled  system  of  the  form 


dt(pi  +  V  •  (v0j)  =  V  •  {oiijV 4>j  +  AiVV20j)  +  giy  i  =  2, 3. 
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where  ay  and  (3j  are  functions  of  0's.  Thus,  the  following  two  equations  are  solved  simulta¬ 
neously, 


30"+1-40?  + 


r1 


2  St 


+  V  •  (v 


n+ 1 


,n+l 


)  =  V- 


a. 


LV0" 


r1) 


+  9 ; 


,n+l 


i  =  2,3. 


The  numerical  scheme  is  implemented  on  a  graphic  processing  unit  using  CUDA.  The 
Nvidia  CUFFT,  Thrust  [03],  as  well  as  CUSP  [□]  have  been  used  to  solve  the  linearized 
system.  The  resultant  solver  is  tested  in  both  time  and  space  to  ensure  it  is  convergent  and 
attains  roughly  second  order. 


4  Numerical  results  and  discussion 

4.1  Parameters  and  initial  setup 

For  convenience,  all  parameters  used  in  the  model  are  listed  in  Table  [0  unless  noticed 
otherwise.  These  parameters  are  chosen  from  the  published  literature  or  user-defined  based 
on  our  best  guesses.  We  note  that,  although  we  have  formulated  the  three-phase  model  by 
treating  the  nucleus  separately  as  a  new  phase  for  eukaryotic  cells,  this  model  is  also  well- 
suited  for  studying  cytoplasmic  dynamics  and  cytokinesis  without  considering  the  nucleus. 
In  those  cases,  we  simply  set  the  initial  condition  of  02  as  zero  and  assume  that  C\  and 
c2  are  the  centers  of  mass  for  the  cytoplasm  distributed  on  each  side  of  the  division  plane, 
respectively.  In  this  two  phase  model,  0i  represents  the  volume  fraction  of  the  substance 
inside  the  cell  membrane  and  03  represents  the  volume  fraction  of  the  buffer,  i.e.  extra 
cellular  matrix  (ECM)  outside  of  the  cell  membrane. 

4.2  Dynamics  of  cell  growth  and  cytokinesis 

During  the  cell  mitotic  process  and  before  cytokinesis  ensues,  the  parent  cell  doubles 
its  cytoplasmic  volume  while  in  the  meantime  duplicates  its  genetic  substances  (DNAs  and 
chromosomes).  The  cell  morphological  change  during  this  process  can  be  studied  using  the 
current  model  by  simply  shutting  down  the  contracting  force  and  switch  on  the  growth 
dynamics  of  the  cell,  which  are  characterized  by  the  moments  known  as  the  checkpoints 
present  in  the  current  model.  At  the  molecular  level,  the  DNA  and  chromosome  replicate 
themselves  and  then  separate  into  two  distinct  sets  of  DNAs  and  chromosomes.  Immediately 
following,  two  offspring  nuclei  form,  each  of  which  contains  the  genetic  information  inherited 
from  the  parent  cell.  This  is  a  complex  process  involving  the  formation  of  mitotic  spindle 
and  its  elongation  in  the  axis  perpendicular  to  the  cytokinetic  ring,  which  is  beyond  the 
scope  of  our  current  multiphase  model’s.  In  the  current  model,  we  focus  on  cytokinesis 
after  the  offspring  nuclei  are  separated.  We  coarse-grain  the  nucleus  as  a  viscous  bulk  fluid 
inside  the  cytoplasm.  After  the  two  offspring  nuclei  are  separated,  the  nuclei  effectively 
position  the  division  plane,  known  as  the  cleavage  plane,  after  which  cytokinesis  ensues. 
When  the  cleavage  plane  is  placed  right  in  the  middle  of  the  long  axis  of  the  cell,  symmetric 
cell  division  can  be  observed.  Otherwise,  cytokinesis  may  result  in  asymmetric  division  or 
failure  of  cytokinesis.  We  will  discuss  the  case  of  asymmetric  division  later,  but  will  not 
discuss  the  failure  of  cytokinesis  since  its  cause  is  not  fully  explored  yet  even  experimentally. 
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In  Figure  □,  we  show  a  detailed  simulation  of  cytokinesis  for  eukaryotic  cells  mentioned 
above.  In  particular,  the  3D  view  of  the  cell  growth  and  division  process  is  portrayed  in 
Figure  [O(a-e).  The  volume  of  the  cytoplasm  is  doubled  at  the  moment  shown  in  Figure  [0(b), 
after  which  we  see  the  formation  of  the  cleavage  furrow  in  Figure  ffl(c),  and  the  abscission 
that  physically  cleaves  the  parent  cell  into  two  offspring  cells  in  Figure  [0(d).  In  Figure  ffl(e), 
we  observe  the  offspring  cell  reshapes  into  round  morphology  under  the  influence  of  surface 
tension  on  its  own  membrane.  To  our  best  knowledge,  this  is  the  first  numerical  simulation 
in  full  3D  using  a  hydrodynamic  phase  field  model.  To  better  visualize  the  dynamics  of 
nucleus  in  the  cytokinesis  process,  a  series  of  2D  slices  at  x  =  0.5  are  plotted  in  Figure  [Q(f-j ) . 
In  particular,  the  cytoplasmic  bridge  connecting  the  two  offspring  cells  is  observed  (shown 
in  Figure  ffl(i),)  which  agrees  quantitatively  with  the  morphogenetic  pattern  obtained  from 
the  experimental  observations  shown  in  Figure  ffl(k). 

We  note  that  the  advantage  of  the  3D  hydrodynamic  model  is  its  capability  to  couple 
the  interior  cytoplasmic  fluid  flow  with  the  exterior  fluid  flow  through  the  cell  membrane 
and  visualize  important  hydrodynamic  quantities  such  as  stress  tensors,  forces  and  the  hy¬ 
drostatic  pressure  throughout  the  domain  during  the  cellular  morphological  transformation. 
In  Figure  □,  the  detail  of  the  viscous  stress  ?/D,  volume-averaged  velocity  v,  as  well  as  the 
hydrostatic  pressure  p  at  time  t  =  8.5  are  depicted.  The  stress  tensor  at  (x,  t)  is  visualized 
as  a  3D  ellipsoidal  object  at  x,  whose  semiaxes  signify  the  length  of  the  three  eigenvalues 
of  the  second  order  symmetric  tensor,  respectively.  From  Figure  E3(b),  a  zoomed  view  of 
viscous  stress  tensor  is  shown.  The  stress  is  highly  inhomogeneous  at  the  contractile  ring 
on  the  division  plane  (cleavage  plane),  which  correlates  well  with  the  velocity  field  shown  in 
Figure  El(f-g) ,  as  the  cytokinetic  ring  of  the  cell  is  contracting.  Besides,  the  gradient  of  the 
hydrostatic  pressure  is  also  high  at  the  cleavage  plane  as  well  as  the  interface  between  the 
cell  nucleus  and  cytoplasm. 

In  addition,  the  distribution  of  the  proxy  force  and  the  surface  tension  force  are  also 
shown  in  Figure  0.  Seen  from  Figure  0(a),  the  proxy  force  is  mainly  distributed  on  the  cell 
membrane  (the  interface  between  cytoplasm  and  ECM),  with  much  higher  values  on  the 
division  plane,  contracting  the  cell  membrane  towards  the  center  of  the  long  axis  of  the  cell. 
This  is  shown  in  the  2D  view  of  proxy  force  in  Figure  0(b-c).  Same  as  the  proxy  force,  the 
surface  tension  force  is  also  distributed  on  the  interfaces.  However,  surface  tension  is  more 
evenly  distributed,  aiming  to  smooth  out  the  interface.  These  two  forces  oppose  to  each 
other  on  the  interface.  It  is  the  competition  between  the  surface  tension  force  and  the  proxy 
force  on  the  division  plane  that  ultimately  contributes  to  the  success  of  cytokinesis. 

4.3  Asymmetric  Cell  Cytokinesis 

The  positioning  of  the  contractile  ring  (or  cytokinetic  ring)  located  on  the  division  plane  is 
dictated  by  the  mitotic  spindle  and  other  intercellular  processes.  It  can  affect  the  morphology 
of  the  offspring  cells  significantly.  If  the  contractile  ring  is  positioned  not  in  a  symmetric 
manner,  then  asymmetric  cell  division,  or  cell  polarization  could  be  observed.  Asymmetric 
cell  division,  which  includes  cell  polarization  and  cytokinesis  is  essential  for  generating  cell 
diversity  during  development. 

Here,  we  conducted  two  case  studies,  where  we  only  concentrate  on  the  cell  membrane 
dynamics  by  setting  02  =  0.  Note  in  this  case,  ci  and  C2  in  the  model  are  chosen  as  mass 
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centers  of  the  half  cell  divided  by  the  cleavage  plane.  For  the  first  numerical  study,  we 
choose  the  contractile  ring  on  the  plane  at  x  =  0.9  with  Lx  =  2.0  for  the  whole  domain.  An 
asymmetric  cell  division  is  simulated  shown  in  Figure  0.  The  asymmetric  cell  division  in  the 
intermediate  stage  resembles  qualitatively  to  the  budding  yeast  polarization. 

Another  case  of  asymmetric  cell  division  is  due  to  the  fact  that  the  concentration  of  actin- 
myosin  filament  is  distributed  inhomogeneously  during  cytokinesis.  As  a  simple  experiment, 
we  set  7^1,  the  contractile  strength  as  a  function  of  space,  to  mimic  this  inhomogenous  con¬ 
tractile  force  due  to  the  heterogenous  distribution  of  actin-myosin  filaments  on  the  division 
plane.  In  this  context,  we  set  7^1  =  7 diy,  be.  the  strength  of  the  contractile  force  is  higher 
on  upper  membrane  of  the  cell.  As  a  result,  a  plane-cell  type  asymmetric  cell  division  can 
be  observed,  shown  in  Figure  0.  The  cleave  furrow  is  observed  in  Figure  0(h),  which  agrees 
qualitatively  with  the  cleavage  in  the  jellyfish  aequorea. 

5  Conclusion 

In  this  paper,  we  use  a  multi-phase  field  model  to  study  cytokinesis  of  an  eukaryotic 
cell  during  its  mitotic  process.  Several  interesting  phenomena  such  as  dynamics  and  mor¬ 
phological  patterns  of  symmetric  or  asymmetric  cell  division  in  cytokinesis  are  numerically 
simulated  with  the  model.  These  morphological  patterns  agree  qualitatively  with  experimen¬ 
tal  observations.  This  simplified  model  is  thus  proven  to  be  an  effective  tool  for  studying 
cytokinesis  during  the  cell  division  process. 

Evidently,  this  simplified  model  needs  significant  improvement  in  order  to  make  it  capa¬ 
ble  to  simulate  the  complex  biological  and  chemical  processes  in  addition  to  the  mechanical 
process  during  the  real  cell  mitotic  process.  Nevertheless,  it  marks  our  Erst  attempt  in  devel¬ 
oping  a  full  3D  hydrodynamic  model  for  cell  division  in  a  hydrodynamically  consistent  way. 
Within  the  framework  for  modeling  cell  dynamics  that  the  model  has  provided,  we  can  fur¬ 
ther  superimpose  additional  features  to  it  or  add  additional  cellular  components/structures 
to  it.  In  our  future  work,  we  will  incorporate  the  actomyosin  microstructure  to  the  cell 
cortical  layer  and  derive  the  contractile  force  on  the  cytokinetic  ring  based  on  their  density 
and  orientational  distribution  instead  of  using  the  proxy  force.  In  addition,  the  viscoelastic 
properties  of  cytoplasm,  the  chromosome  spindle,  as  well  as  its  elongation  on  facilitating 
cytokinesis  could  be  added  incrementally  as  well. 
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Tabic  1:  Dimensional  Parameters. 


Symbol 

Description 

value 

Unit 

Reference 

T 

Absolute  Temperature 

303 

Kelvin 

PH 

k 

Boltzmann  constant 

1.38  x  10~23 

m2kgs~2  K-1 

PH 

h 

Characteristic  length  scale 

5  x  10~5 

m 

ra 

to 

Characteristic  time  scale 

0.1 

s 

m 

p 

Reference  density 

1.1  x  103 

kgm~ 3 

PH 

m 

Dynamic  viscosity  of  cytoplasm 

1  x  10~2 

kgm~  Ls~ 1 

m 

V2 

Dynamics  viscosity  of  nueclues 

2  x  10-2 

kgm~  1s~ 1 

m 

Dynamics  viscosity  of  ECM 

5  x  10-3 

kgm~  1s~ 1 

m 

1A,B,C 

Distortional  energy  coefficient 

2.5  x  10Y 

m~ 1 

71,2,3 

Bulk  free  energy  coefficient 

1.875  x  102° 

m~3 

7123 

Lagrangian  multiplier  strength 

1.0 

m~ 3 

Al,2,3 

Motility  parameter  for  each  component 

1  x  10"11 

kg~1m^s 

Cl 

Growth  rate  of  cytoplasm 

5.0 

S~L 

C2 

Growth  rate  for  nucleus 

5.0 

s-1 

7dl 

Stimulating  force  strength 

2.5 

m~ 3 

7d2 

Stimulating  force  strength 

15 

m~3 

Sd 

0.01 

15 


(f)  t  =  0  (g)  t  =  1.25  (h)  t  =  2.5 


(i)  t  =  15  (j)  t  =  26.5  (k) 


Figure  1:  Cell  growth  and  cytokinesis.  This  figure  shows  a  parent  cell  duplicates  its 
nucleus  and  cytoplasm,  then  splits  into  two  identical  offspring  cells,  (a)-(e)  3D  Nu¬ 
merical  simulations  of  the  cell  division  process  at  different  stages;  (f-j)  2D  slices  at 
x  =  0.5  for  the  cell  division  process;  (k)  a  dividing  melanoma  cell  just  before  it 
divides  into  two  offspring  cells  completely.  Permission  is  needed  from,  http  : 
// php.med.unsw .edu.au / cellbiology / index. phpltitle  =  File  :  Melanomacytokinesis.jpg). 
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(a)  ?/D 


(b)  Zoom  view  of  77D 


(c)  77D  at  x  =  0.75  (d)  Zoom  view  of  77D  at  x  =  0.75 


(e)  77D  at  z  =  0.5 


I 


0.050 

0.038 

0.025 

>- 

0.013 

7.9e-18 
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(f)  vati  =  0.75 


(g)  v  at  x  =  0.5 


(h)  p  at  x  =  0.75 


Figure  2:  Hydrodynamic  variables  at  t  —  8.5.  This  figure  shows  the  hydrodynamic  variables, 
including  viscous  stress  tensor  (3 D,  volume-averaged  velocity  v  and  the  hydrostatic  pressure 
p  distribution  at  time  t  =  8.5  for  the  simulation  shown  in  Figure  □.  (a)  3D  view  of  the  stress 
tensor  as  ellipsoids;  (b)  a  zoomed  view  of  the  stress  tensor;  (c-e)  2D  slices  of  the  stress  tensor 
field;  (f-g)  2D  slices  of  the  velocity  filed;  (h)  a  2D  slice  of  the  pressure  field. 
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Figure  3:  Visualization  of  the  proxy  force  and  the  surface  tension  force.  This  figure  shows 
the  proxy  force  and  surface  tension  force  at  time  t  =  8.5  of  the  simulation  shown  in  Figure 
□.(a)  3D  view  of  the  proxy  force;  (b)  a  2D  slice  of  proxy  force  at  x  =  0.5;  (c)  a  2D  slice  of 
proxy  force  at  z  =  0.5;  (d)  3D  view  of  the  surface  tension  force;  (e)  a  2D  slice  of  the  surface 
tension  force  at  x  =  0.5;  (f)  a  2D  slice  of  the  surface  tension  force  at  z  =  0.5.  The  proxy 
force  opposes  to  the  surface  tension  force  at  the  cytokinetic  ring. 
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(a)  t  =  0 


(b)  t  =  5 


(e)  t  =  0  (f)  t  =  5  (g)  t  =  10 


oO 


(h)  t  =  18 

Figure  4:  Asymmetric  cell  division  due  to  the  asymmetric  positioning  of  cleavage  plane.  This 
figure  shows  an  asymmetric  cell  cytokinesis  process  in  which  the  division  plane  is  positioned 
in  an  asymmetric  fashion  along  the  long  axis  of  the  cell,  (a)-(d)  3D  Numerically  simulated  cell 
division  process  at  different  time;  (e-h)  2D  slices  ( z  =  0.5  )  of  cell  division  process  at  different 
time;  (f)a  budding  yeast.  This  asymmetric  cell  division  resembles  the  yeast  cell  budding 
process.  A  permission  is  needed  from  http:  // www.ppdictionary.com/ fungi.htm 
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(a)  t  =  0 


(c)  i  =  15 


(b)  t  =  10 


(d)  t  =  20 
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(e)  t  =  0 


(f)  t  =  10 


>4 


(g)  t  =  15  (h) 

Figure  5:  Asymmetric  cleavage  furrow  formation  due  to  inhomogeneous  actomyosin  dis¬ 
tribution.  This  figure  shows  the  asymmetric  cleavage  furrow  formation  due  to  the  in- 
homogenous  contractile  force  along  the  cytokinetic  ring  induced  by  heterogenous  acto¬ 
myosin  distribution  on  the  contractile  plane,  (a-d)  3D  view  of  the  asymmetric  cell  divi¬ 
sion  at  different  time;  (e-g)  2D  slices  (z  =  0.5)  of  the  asymmetric  cell  division  at  differ¬ 
ent  time;  (h)  cleavage  furrow  in  the  jellyfish  Aequorea.  A  permission  is  needed  from 
http  :  I/ celldynamics.org/embryos/aequore2Qitml) . 


